C********************************************************************

      SUBROUTINE USRCRP(MODEL,FF,DFDS,EFSTS,TEMP,TEPS)

C********************************************************************
C
C     THIS SUBROUTINE CALCULATES THE CREEP RATE AND ITS DERIVATIVE.
C
C     INPUT :
C
C      MODEL     =   CREEP EQUATION NUMBER
C      EFSTS     =   EFFECTIVE STRESS
C      TEMP      =   TEMPERATURE
C      TEPS      =   EFFECTIVE STRAIN (ONLY USED TO EVALUATE THE 
C                    CREEP RATE, NOT TO AFFECT THE FLOW STRESS)
C
C     OUTPUT :
C
C      FF        =   CREEP RATE
C      DFDS      =   DERIVATIVE OF CREEP RATE W.R.T. STRESS
C
C********************************************************************

      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C  USER SUPPLIED SUBROUTINES
C
C ****    USER DEFINED VARIABLES ****
C
      CHARACTER*80 IUSRVL
      COMMON /IUSR/ IUSRVL(10)
C     
C     TO READ DATA (10 RESERVED LINES)
C        READ(IUSRVL(LINE NUMBER),*) DATA1,DATA2,DATA3...
C
C     TO WRITE DATA (10 RESERVED LINES)
C        WRITE(IUSRVL(LINE NUMBER),*) NEWDATA1, NEWDATA2, NEWDATA3 ...
C
C ****     END    ****
C
C     Version 6.0
C
C     COMMON /ELMCOM/
C
C         RZE     : Four corner coordinates
C         URZE    : Velocity
C         STSE    : Stress
C         EPSE    : Strain rate
C         EFEPSE  : effective strain rate
C         EFSTSE  : Effective stress
C         TEPSE   : Total effective strain
C         RDTYE   : Density
C         TEMPE   : Temperature
C         DTMPE   : Temperature rate
C         DAMAGE  : Damage value
C         USRE1   : Element user variable array at the begining of step
C         USRE2   : Element user variable array at the end of step
C         USRNE   : Nodal user state variables at 4 nodes
C         NODEE   : Connectivity
C         KELE    : Global element number
C         KELEL   : Local element number
C         KGROUP  : Material group number
C
      COMMON /ELMCOM/ RZE(2,4),URZE(2,4),STSE(6),EPSE(6),EFEPSE,EFSTSE,
     +                TEPSE,RDTYE,TEMPE(4),DTMPE(4),DAMAGE,
     +                USRE1(1500),USRE2(1500),
     +                USRNE(1500,4),NODEE(4),KELE,KELEL,KGROUP

C
C     COMMON /NODCOM/
C
C        RZN      : Nodal point coordinates
C        URZN     : Nodal point velocities
C        DRZN     : Nodal point displacement
C        TEMPN    : Nodal point temperature
C        USRN1    : User defined state variables (Input : At the beginning of Step N)
C        USRN2    : User defined state variables (Output: At the end of Step N)
C        KNODE    : Node number
C
      COMMON /NODCOM/ RZN(2),URZN(2),DRZN(2),TEMPN,DTMPN,USRN1(1500),
     +                USRN2(1500),KNODE

C
C     CURTIM: CURRENT TIME           
C      
      COMMON /CLOK/ CURTIM
C
C     DTMAXC: CURRENT TIME STEP SIZE
C
      COMMON /SSTU/ DTMAXC
C
C
C     COMMON /USRCTL/
C
C        KOBJ      : Object number
C        KSTEP     : Step number (N)
C        ISTATUS   : 0 - the beginning of the step
C                    1 - the end of the step
C
C     WHEN (ISTATUS.EQ. 1)  --> USRE2/USRN2 should be updated here
C     KELE  > 0             --> Element data is active
C     INODE > 0             --> Node Data is active
C
      COMMON /USRCTL/ KOBJ,ISTATUS,KSTEP,KSSTEP
C
C     Branching to proper creep rate routine based on the
C     number specified in the pre-processor
C     
      GO TO (510,520,530,540,550), MODEL 
C     
  510 CALL UCRP1(FF,DFDS,EFSTS,TEMP,TEPS)
      RETURN
C
  520 CALL UCRP2(FF,DFDS,EFSTS,TEMP,TEPS)
      RETURN
C
  530 CALL UCRP3(FF,DFDS,EFSTS,TEMP,TEPS)
      RETURN
C
  540 CALL UCRP4(FF,DFDS,EFSTS,TEMP,TEPS)
      RETURN
C
  550 CALL UCRP5(FF,DFDS,EFSTS,TEMP,TEPS)
      RETURN
C
C  TO BE CONTINUED BY THE USER, IF NECESSARY
C
      END
C
C********************************************************************
C
      SUBROUTINE UCRP1(FF,DFDS,EFSTS,TEMP,TEPS)
C
C********************************************************************
C
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C     DESCRIPTION
C
C       THIS ROUTINE IS USED TO DEMONSTRATE THE IMPLEMENTATION OF THE
C       CREEP RATE ROUTINE.  ALL THE REAL VARIABLES SHOULD BE IN DOUBLE 
C       PRECISION. THE DEFINITION OF THE ABOVE PARAMETERS ARE DESCRIBED 
C       AS FOLLOWS.
C     
C     INPUT :
C
C      MODEL     =   CREEP EQUATION NUMBER
C      EFSTS     =   EFFECTIVE STRESS
C      TEMP      =   TEMPERATURE
C      TEPS      =   EFFECTIVE STRAIN (ONLY USED TO EVALUATE THE 
C                    CREEP RATE, NOT TO AFFECT THE FLOW STRESS)
C      VALUES(NGRP,NTMP)
C                =   TEMPERATURES AND TEMPERATURE-DEPENDENT VARIABLES
C                    FOR THE CREEP RATE EQUATION
C      NTMP      =   NUMBER OF TEMPERATURES
C      NGRP      =   NUMBER OF TEMPERATURE-DEPENDENT VARIABLES + 1         
C
C     OUTPUT :
C
C      FF        =   CREEP RATE
C      DFDS      =   DERIVATIVE OF CREEP RATE W.R.T. STRESS
C
C     CREEP RATE MODEL USED IN THIS EXAMPLE:
C
C      FF        =   GAMMA * (EFSTS) ^ EM
C      DFDS      =   GAMMA * EM * (EFSTS) ^ (EM - 1)
C      
C        WHERE "GAMMA" AND "EM" ARE TEMPERATURE-DEPENDENT AND ARE
C        GIVEN IN A TABLE FORM.  THEIR ACTUAL VALUES ARE TO BE
C        INTERPOLATED BASED ON THE CURRENT TEMPERATURE "TEMP".
C
      DIMENSION VALUES(3,9),INDX(2),RES(2)
C
      IF(EFSTS.EQ.0.) THEN
        FF = 0.
        DFDS = 0.
        RETURN
      ENDIF
C
      NTMP = 9
      NGRP = 2 + 1
C
C  ENTER TEMPERATURES IN THE FIRST ROW OF "VALUES"
C
      VALUES(1,1) = 10.0e0
      VALUES(1,2) = 200.0e0
      VALUES(1,3) = 400.0e0
      VALUES(1,4) = 600.0e0
      VALUES(1,5) = 800.0e0
      VALUES(1,6) = 1000.0e0
      VALUES(1,7) = 1200.0e0
      VALUES(1,8) = 1400.0e0
      VALUES(1,9) = 1600.0e0
C
C  ENTER "GAMMA" IN THE SECOND ROW OF "VALUES"
C
      VALUES(2,1) = 1.14e-13
      VALUES(2,2) = 1.14e-13 
      VALUES(2,3) = 1.14e-13
      VALUES(2,4) = 1.14e-13
      VALUES(2,5) = 1.14e-13
      VALUES(2,6) = 1.14e-13
      VALUES(2,7) = 1.14e-13
      VALUES(2,8) = 1.14e-13
      VALUES(2,9) = 1.14e-13
C
C  ENTER "EM" IN THE THIRD ROW OF "VALUES"
C
      VALUES(3,1) = 4.80e0
      VALUES(3,2) = 4.80e0
      VALUES(3,3) = 4.80e0
      VALUES(3,4) = 4.80e0
      VALUES(3,5) = 4.80e0
      VALUES(3,6) = 4.80e0   
      VALUES(3,7) = 4.80e0   
      VALUES(3,8) = 4.80e0   
      VALUES(3,9) = 4.80e0   
C
C  IN "INDX(1:2)" ENTER THE METHOD OF INTERPOLATION: 
C       = 0 -- LINEAR INTERPOLATION
C       = 1 -- LOGARITHMIC INTERPOLATION
C
      INDX(1) = 1
      INDX(2) = 0
C
C  CALL "INTPOL4", A "DEFORM" INTERPOLATION SUBROUTINE 
C
      CALL INTPOL4(RES,TEMP,VALUES,NTMP,NGRP,INDX)
C
C  FROM RES(1:2) THE CONSTANTS ARE OBTAINED BASED ON THE CURRENT TEMPERATURE "TEMP"
C
      GAMMA = RES(1)
      EM = RES(2)
C
C  CALCULATE "FF" AND "DFDS"
C
      FF = GAMMA * EFSTS**EM
C      DFDS = GAMMA * EM * EFSTS**(EM - 1.)
      DFDS = FF * EM / EFSTS
C
      RETURN
      END
C********************************************************************
C
      SUBROUTINE UCRP2(FF,DFDS,EFSTS,TEMP,TEPS)
C
C********************************************************************
C
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C     SEE UCRP1 FOR DETAILED DESCRIPTION
C
      RETURN
      END
C********************************************************************
C
      SUBROUTINE UCRP3(FF,DFDS,EFSTS,TEMP,TEPS)
C
C********************************************************************
C
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C     SEE UCRP1 FOR DETAILED DESCRIPTION
C
      RETURN
      END
C********************************************************************
C
      SUBROUTINE UCRP4(FF,DFDS,EFSTS,TEMP,TEPS)
C
C********************************************************************
C
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C     SEE UCRP1 FOR DETAILED DESCRIPTION
C
      RETURN
      END
C********************************************************************
C
      SUBROUTINE UCRP5(FF,DFDS,EFSTS,TEMP,TEPS)
C
C********************************************************************
C
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
C
C     SEE UCRP1 FOR DETAILED DESCRIPTION
C
      RETURN
      END

C **********************************************************************
      SUBROUTINE INTPOL4(RES,XX,VAL,NX,ND,INDX)
      IMPLICIT REAL*8 (A-H,O-Z), INTEGER*4 (I-N)
      DIMENSION VAL(ND,*),RES(*),INDX(*)
C
C  THIS IS THE INTERPOLATION ROUTINE FOR ONE VARIABLE AND 
C  THE DATA ARE STORED IN A MATRIX OF (ND,NX)
C
C  INDX(I) = 0: LINEAR INTERPOLATION
C          = 1: LOGARITHMIC
C
      NGRP = ND - 1
      IF(XX.LE.VAL(1,1)) THEN
        DO 20 I = 1, NGRP
        J = I + 1
   20   RES(I) = VAL(J,1)
        GOTO 100
      ENDIF
C
      DO 30 I = 2, NX
      IF(XX.LE.VAL(1,I)) THEN
        NSEGA = I - 1
        NSEGB = I 
        GOTO 50
      ENDIF
   30 CONTINUE
      NSEGA = NX - 1
      NSEGB = NX
C
   50 X1 = VAL(1,NSEGA)
      X2 = VAL(1,NSEGB)
      DL1 = XX - X1
      DLA = 1.D0/(X2 - X1 + 1.D-20)
      R = 2.0D0*DL1*DLA - 1.0D0
      H1 = 0.5D0*(1.0D0 - R)
      H2 = 0.5D0*(1.0D0 + R)
C
      DO 70 I = 1, NGRP
      J = I + 1
      F1 = VAL(J,NSEGA)
      F2 = VAL(J,NSEGB)
      IF(INDX(I).EQ.0) THEN
        RES(I) = H1*F1 + H2*F2
      ELSE
        IF(F1.GT.0.) THEN
          FLG1 = LOG(F1)
        ELSE
          FLG1 = -1.D10
        ENDIF
        IF(F2.GT.0.) THEN
          FLG2 = LOG(F2)
        ELSE
          FLG2 = -1.D10
        ENDIF
        DUM = H1*FLG1 + H2*FLG2
        RES(I) = EXP(DUM)
      ENDIF
 70   CONTINUE
C
  100 RETURN
      END
C---+----*----+----*----+----*----+----*----+----*----+----*----+----*--
